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Abstract 



The determination of a steeply falling energy spectrum from rare events with non- 
gaussian measurement errors is a delicate matter. The final shape of the spectrum may 
^ . be severely distorted as a consequence of non-gaussian tails in the energy resolution 

fin \ of the experiment. We illustrate this effect with the recent experimental efforts to 

qq ' determine the ultra-high energy extreme of the cosmic ray spectrum. 
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^ ■ 1 Introduction 

We will study the determination of a steeply falling spectrum when the energy measurement 
errors show long, asymmetrical, non-gaussian tails. We will see that for a steep spectrum the 
effect of long tail error functions can be dramatic: spectrum slopes, cutoff or other features 
could be changed or even erased, depending on the degree of deviation of the error function 
(energy resolution) from a gaussian distribution. 

This problem, well known to nuclear spectroscopists working near the end point of a 
spectrum |l| is, however, often ignored in areas such as Cosmic Ray physics and we there- 
fore choose to illustrate it by examining the spectrum of Ultra High Energy Cosmic Rays 
(UHECR) as determined by ground arrays, addressing the crucial question of the existence 



or absence of a cutoff in the spectrum, known as the Greisen-Zatsepin-Kuzmin (GZK)cutoff 



2 Spectrum reconstruction with non-gaussian errors 

We will consider four cases for the reconstruction error function. These cases illustrate the 
different possibilities. 

Let's start with the simplest, well-known, case of a gaussian error reconstruction function. 
Assume that events of energy E' are reconstructed with energy E with a probability given 
by 

P(E,E>) = NeM ~ iE ~ E ' )2 ), (1) 

where a is constant. Then, the effect on an initial power law spectrum, = AE" 1 , will be 
simply given by the convolution of <fi with the error function 

roc 

<p{E) = dE'P(E,E')4>(E'). (2) 
Jo 
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In this case the integral can be carried out analytically but is not very illuminating. For 
high energy E it can be expanded to give 



4>(E) =A2T^1 + 7(7 + ...)• (3) 

At large energies, the effect of the reconstruction energy on the reconstructed flux is very 
small. This is what is expected for gaussian errors and therefore one can safely neglect this 
effect. The same occurs for any other initial spectrum. The reconstructed spectrum will 
differ from the original one by a function which goes to zero rapidly with increasing energy. 
Now let's consider the log-gaussian distribution 

mm = N^J^l), (4, 

where N is a normalization constant and A£ is the standard deviation. Then the convolution 
of P with the power like initial flux is 

4>(E) = [°° dE'<j>{E')P{E,E') = A£^exp(^( 7 2 - 1)), (5) 
Jo 2 

where we are assuming that A£ is independent of energy. We see therefore that for a log- 
gaussian error with a constant standard deviation the effect of the convolution is to change 
the normalization. This is rather different from the result of Eq.^|, where the change was 
negligible at large enough energies. Here the flux increases (for 7 > 1) by a constant factor. 
Notice also that the enhancement factor depends strongly on the initial spectrum index, 7. If 
we put 7 = 3 and A£ ~ 15% then 0/0 ~ 1.09, for A£ ~ 30% then 0/0 ~ 1.43. Alternatively, 
for 7 = 3, we get a factor two enhancement for A£ ~ 0.4. 

In order to see what is the effect of very long tails on the reconstruction of events, let us 
consider an error function with a long tail 



P(E,E') = N(^r(l + ^)- a -^ (6) 



where N is a normalization and a and ft are constants. The energy is reconstructed with a 
power law distribution of slope a and —ft for small and large reconstructed energies. The 
convolution with the original spectrum gives as before 

4 > (E)=AK(a,(3, 1 )E-\ (7) 

where K is a constant which depends on the indexes a, ft, 7. As in the log-gaussian case the 
normalization changes but now the effect can be dramatic, depending on the values of a and 
ft. For a = 4 and ft = 7 we find 0/0 ~ 1.5. The log-gaussian distribution can be thought as 
the limit when a and ft are large. 

In both cases the normalization is strongly affected. This is easily understood, the long 
tails on the error functions favor the more numerous low energy events to be reconstructed 
at higher energy, giving a higher flux. In all these cases we are assuming that there are 
no systematical energy shifts: the average reconstructed energy for a fixed initial energy 
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is assumed to be equal to the initial energy. The existence of systematic errors in the 
reconstruction of energy would worsen considerably our results since in this case the change 
on the normalization of the spectrum will be linear in A£ syst rather than quadratic. Also we 
are assuming that the energy reconstruction function has no strong energy dependence. 

Our last example of an error function with long tails is the Moyal distribution which 
constitutes an approximation to the Landau distribution describing the fluctuations in energy 
losses of an ionizing particle passing through matter 

1 (E — E') (E-E>) 

P(E, E') = Nexp(--( { - >- + e-M . (8) 

I a 

We will return to this distribution at the end of the paper when we look in detail at the 
problem of the GZK cutoff in the UHECR spectrum. Now let's consider an initial spectrum 
with some feature. In order to keep the analysis simple, let's assume a spectrum with a 
exponential cutoff. 

^E) = AE^exp(-^), (9) 

where E c is the cutoff energy. The effect of a gaussian error would be again the original 
flux times a function going rapidly to 1 with energy, i.e. the original flux is not modified. 
However for a log-gaussian error or for a power law error we will have strong modification 
of the flux. In both cases at low energy (E <C E c ) we will have the previous result, the 
reconstructed flux is modified by a constant factor. At high energies, E E c , the effect is 
different. In the log-gaussian error function we have 

4>(E) = AE' 1 exp(A^/2( 7 2 - l))e" B/ ^, (10) 

where now E c = i? c exp(7A£ 2 ) is the new, smeared, cutoff energy. Again, this effect can be 
large due to the non linear dependence on 7 and A£; for A£ = 0.3 and 7 = 3 E c = 1.3E C , 
for A£ = 0.4 E c = l.62E c . 

In the case of a power law error function the effect is even more dramatic. At high energy, 
E ^> E c , the convolution of the error function with the exponential cutoff flux is 

4>(E) = AE^Et"K\a } f3 }7 ), (11) 

i.e. the cutoff is completely washed out and replaced by a spectrum with a constant slope. 
The new slope is solely dependent on the error function and not on the original spectrum. 

We want now to argue that the non-gaussian error functions given above are realistic. 
There are several reasons why one expects long tails in the reconstruction of events in cosmic 
ray experiments. In ground array experiments the energy is estimated by measuring the 
density of particles in the shower. Usually parameters such as p(600), the density of particles 
at 600 m from the core, are used. It is found that the energy of the shower scales as 
E oc p(600) K and k ~ 1. But this quantity has non-gaussian uncertainties. First, to calculate 
p(600) (or any similar parameter) the core position and the arrival direction of the shower are 
needed. Given an arrival direction the p(600) is corrected to "zero degrees" arrival direction 
by an equation like 

p(600) = p(600) exp(s (sec(fl) - 1)), (12) 
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Figure 1: Distribution of the total number of charged particles at ground level for protons 
of E = 1 EeV generated using Aires with Sibyll hadronic generator. 
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Figure 2: Distribution of p(600) at ground for protons of E = 1 EeV generated using Aires 
with Sibyll hadronic generator. 
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where 6 is the zenith angle and So is a constant. Even if the arrival direction is reconstructed 
with a gaussian probability function the resulting p(600) is not. However, this effect is, 
we believe, small at all but the largest zenith angles. The effect on the core reconstruction 
is more important. The density of particles at a distance r from the shower core can be 
parameterized by the NKG lateral distribution function 

p(r) = Kr\l + r/Ry', (13) 

where K, R, rj, rj' are constants. Then, the error in the density due to the error in the 
determination of the core position would be generally non-gaussian. 

Finally, any indirect measurement of cosmic rays based on shower development is sub- 
jected to shower fluctuations. It is well known that the total number of particles at fixed 
depth have large non-gaussian fluctuations, the total number of particles r.m.s scales as the 
energy. The same occurs for p(600), since it is related to the total number of particles. In 
figure p] we can see the distribution of the total number of particles at ground level in a 
shower of fixed energy. We can appreciate that the distribution of the number of particles is 
far from a gaussian and has long tails. In figure H we show the distribution of the density at 
600 meters for showers of fixed energy. More than 20000 protons showers of energy 10 18 were 
simulated with the Aires Monte Carlo code |4[]. As before the distribution is non-gaussian 
with tails extending to more than 1 order of magnitude above the average. 

This was known for a long time but rarely used in the analysis of cosmic ray physics 
||. Apparently, it is a generic phenomena in probability theory and does not contradict the 
central limit theorem 0. 

3 The GZK Cutoff 

Our results are not directly applicable to any real experiment. Rather, experiments should 
ascertain what are their error reconstruction functions and take them into account in the 
calculation of the spectrum. Also, events with large fluctuations could be cut off by other 
methods, improving therefore the energy reconstruction and avoiding these unwanted effects. 
Particularly important would be to get rid of any power like tail in the reconstruction error 
functions, since they have the most devastating effects on the spectrum reconstruction. On 
the other hand this effect could, at least in part, explain the current experimental situation 
where the normalization of the cosmic ray flux for different experiments is different and the 
existence of the GZK cutoff is controversial. Values of A£ ~ 30 — 40 % are usual in cosmic 
ray physics and as said before with such values there are important effects on the shape and 
normalization of the resulting spectra. 



5 



I . . ..I 



20 50 
E r (EeV) 



100 



200 



Figure 3: Expected flux for log-gaussian energy reconstruction error function with A£ = 
0.1,0.3,0.6 (continuous, dashed, and dotted lines) for a flux with exponential cut-off at 
E = 50 EeV and for a cosmological distribution of sources. 

The current experimental situation is aggravated by the low statistics in the high energy 
part of the spectrum. The Agasa experiment has measured 58 || events above 4 x 10 19 eV, of 
them 8 are above 10 20 eV. Using a realistic spectrum obtained from an uniform cosmological 
distribution of sources || and using a log-gaussian error function with A£ = 0.4 we obtain 
that for a total of 58 events above 4 x 10 19 eV, the probability of having 8 or more events 
above 10 20 eV is 2%. The average number of expected events would be 3. In figure ^ the 
convolution of different spectrum with energy error functions is shown. In the case of a 
A£ = 0.3 and the cosmological distribution of sources of Teshima and Yoshida || one can 
see that at E — 200 EeV only a 0.1 reduction on the flux is expected. 

The number of AGASA events above the GZK cutoff seems to establish its absence. 
However, if we introduce non-gaussian errors, we stress once again, the observation of such 
events becomes compatible with the GZK cutoff. To illustrate quantitatively the possibility 
of not seeing the GZK cutoff because of non-gaussian errors in the energy determination we 
consider the Moyal distribution of equation 4 (with approximately the same half-width of a 
gaussian distribution of a given variance, as seen in fig. 4) to each energy sampled from the 
AGASA spectrum with a exponential cutoff, that is a power law spectrum with an index of 
2.78 multiplied by an exponentially falling spectrum starting at 4 x 10 19 eV. These errors 
alter the shape of the original spectrum and the expected number of events above the GZK 
cutoff . We summarize the results of this analysis in fig. 5. We can compare the probability 
of finding events above GZK for measurements affected by gaussian errors (10%, 20%, and 
30%) with measurements affected by Moyal-distributed errors (with the same half-width of 
gaussians with standard deviation of 10%, 20%, and 30% of energy). We see that considering 
only gaussian errors we cannot explain the observed number of events above GZK by AGASA 
(indicated by an arrow in figure 5) and would be led to conclude that there is a cutoff in 
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Figure 4: The Moyal distribution of Eq. 4 compared to a gaussian distribution having 
approximately the same half width 
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Figure 5: Probability of obtaining events above the GZK cutoff for different resolutions (error 
functions) 
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the UHECR spectrum. However, for Moyal-distributed errors (30% of energy) we have an 
appreciable probability of finding more than 47 events above the GZK cutoff. We have 
shown in this paper that the experimental resolution curve when deviating from a gaussian 
distribution severely impacts the shape and normalization of a steeply falling spectra near its 
end point. Spectral features such as a cutoff may disappear from the reconstructed spectrum 
due to the effect of long tails in the resolution(error)function of the experiment. 
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